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COMMENTS ON THE ROLE OF DIAGONAL DOMINANCE 


IN IMPLICIT DIFFERENCE METHODS 


By David H. Rudy euid Richard S. 

Langlej Research Center 


Hirsh 


SUMMARY 


Numerical tests were made for a model of the Navier-Stokes equations 
using a second-order accurate implicit scheme i^ich guarantees diagonal 
dominance. The results of these tests suggest that the failure of implicit 
methods using large marching steps may not always be attributed to the lack of 
diagonal dominance in the coefficient matrix. In some cases the failure may 
be caused by a nonlinear instability associated with the solution method. 

INTRODUCTION 


Implicit finite difference methods in genersJ. offer the major advantage of 
unconditional stability (with respect to the step size in the marching direc- 
tion) over explicit methods. This advantage is offset somewhat by the need to 
solve a system of simultaneous algebraic equations during each marching step. 
However, when central differences are used, the coefficient matrix becomes 
tridiagonal in form and can be inverted much more rapidly than a full matrix. 
Consequently, to gain a computational time advantage, implicit marching steps 
need only to be a factor of four or five times larger than the maximum step 
size allowed in an explicit method. In many cases, though, large marching 
steps cause the coefficient matrix to lose diagoneil dominance, a sufficient, 
but not necessary, condition for convergence of the matrix inversion. Without 
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diagonal dominance, roundoff error can accumulate in the matrix inversion and 
destroy the solution. For Burgers* equation (a model of the Navier-Stokes 
equations), Hirsh and Rudy (Ref. l) fo\uid that diagonal dominance could he 
maintained for any size marching steps if the cell Reynolds number was 2 or 
less. In the present paper results are given for numerical tests made with 
the second-order acc\irate scheme of Khosla-Rubin (Ref. 2) which maintains 
diagoncd dominance for all cell Reynolds numbers and meui'ching step sizes. 

The results of these tests suggest that the failure of implicit methods using 
large marching steps may not always be caused by the lack of diagonal dominance. 


a 

A,B,C,D 

C 

D* 

R 

c 

u 

U 

Ax 

At 

c 

V 

Subscripts 

i 

max 

x,t 



index denoting grid point spatial location 
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derivative with respect to x-direction, time 
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Superscripts 

n index denoting time level 

DESCRIPTION OF NUMERICAL TESTS 

The numerical tests of the Khosla-Rubin scheme were first made using 
the nonconservative form of Burgers’ equation 


u. + au„ » vu 
t X XX 


( 1 ) 


(where a can he a function of u, t, and x) for several cases presented in 
Ref, 1 where diagonal dominance was lost. The difference fonn of eq. (l) 
becomes, using the Khosla-Rubin scheme (ref. 2), 
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where a^ = u^ - U and U is the steady-state wave speed. Thus, 


* ®i“i * = ’’i 
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The boundary conditions were taken to be: 
u(x,t) = 1.0, at X = -5 
u(x,t) = 0, at X = +5 

Fifty-one equally spaced mesh points were used with Ax = 0.2. The wave speed 
U was 0.5. A linear initial velocity distribution between the boundary 
points was used. The Thomas algorithm was employed to solve the resultant set 
of simultaneous algebraic equations at each time step. The solution was 
assumed to have reached a steady-state when the maximum change in u between 
time steps was less than some specified value t. 

RESULTS AND DISCUSSION 

Nonconservative Form of Nonlinear Burgers' Equation 

For V = 1/24 (R = 2.4) and C « = 10, the solutions converged 

c max 2Ax 

using the Khosla-Rubin scheme in less than 50 steps to the steady-state 
solution presented in Ref. 1 for e = lO”^. The solution is shown in 
Figure 1. In the calculation, the nonlinear coefficient, a^^, was lagged one 
time-step and no iteration was performed at each time step. The solution was 
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then continued further in an attempt to satisfy e * After 

approximately 30 more steps the solution began to diverge rapidly from the 

previous result. (Converged results were obtained for e * and large 

values of C with the usual central difference form.) Similar behavior 
max 

was observed for C * 0.9. For C “0,1, the solution converged in 7^*1 
max max 

steps for e “ lO”^^. In all cases the coefficient matrix is always diagonally 
dominant. All test results given in the paper are summarized in Table I. 

An effect of this instability in the solution can be seen in figure 2 

where 


iDil 


I ^i-t-l ~ "^i ^i-l | 

' 2Ax ' 


is plotted for i “ 26 (the midpoint) as a function of the number of time steps 

for C = 0.1, 0.9 and 10.0. This term, Dl , is the quantity added to the 
max 1 

first-order (upwind) differencing of u^ in the Khosla-Rubin modification and 

appears on the right-hand side of (3), not as a coefficient of the unknowns. 

As shown in figure 2 , D ' . , remains small for C =0.1; however , for 

i=«'D max 

*^max ” 10.0, the term grows rapidly, and large errors appear in u. 

Two types of iteration of the nonlinear term were also investigated for 

the nonconservative form of Burgers' equation. The first of these was the 

use of multiple iterations at each time step in the procedure previously 

described. The case with v = l/2l; and C = 0.9 was run with a fixed 

max 

number of iterations to update a^ at each time step. However, the solutions 
again rapidly diverged with 2, U, or 8 iterations per time step. 

The second form of iteration used was Newton iteration (ref, 3). In this 


technique 
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n+l,m+l n+l,m ^ . 

Ui ’ « 4- 6Uj^ 

where m is the iteration index a'", each time step. Equation (5) is substi- 
tuted into equation (2) to give tridiagonal equations of the form 


♦ Si*"! ♦ 


( 6 ) 


where 


, for u^ < 1/2 
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and for u^ > 1/2, 
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(Ax) 


1 / n+l,m , /o\/ n+l,m n-*-l,n. 
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It should he noted that th-^ convective derivative in has been evaluatc»d 

using a central difference with values at the last iterate of the present time 
step rather than with the usual Khosla-Rubin term. This modification was 
found to be necessary to maintain the asymmetry of tie solution about the 
midpoint. 

Calculations were made in which the solution wa.i iterated at each time 

step until the difference between iterations was less than lO”^*^. (The Newton 

procedure could also be performed with no multiple iteration at each time 

step, but such an approach was not used in the present investigation.) The 

steady-state convergence criterion was again e « lO”^®. Converged solutions 

were obtained for v = 1/21+ and C * 0.9* However, the solution diverged 

max 

when C was increased to 10. 
max 


Linear Burgers' Equation 

To further demonstrate that the observed divergence is related to the 
nonlinear term, tests were also made using a linear Burgers' equation 


u. + au = vu 

t X XX 


( 7 ) 


where a is a constant. The bo’ondary conditions were again taken to be 
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u(x,t) ■ 1.0 at X ■ -5 
u(x,t) *0 at X ■ '♦■5 


For this solution domain, the analytical solution is 


5w 

e - e*^ 

-5« 5w 

e - e 


( 8 ) 


where w » a/v. For v » 1/2U, equation (8) indicates that on extremely steep 
gradient exists which could not be resolved with the coarse Ax » 0.2 grid; 
therefore, v » 1 was used with values of a chosen such that > 2. Steady 
state solutions are shown in figure 3 for two cases, ■ 2.4 and R^ ■ 10. 
Convergence was obtained in both cases for e w 10~^^ for values of C up 
to 1000. The "wiggles” in the solution are not roundoff error but enre the 
exact algebraic solution of the difference equations. 

Conservative Form of Nonlinear Burgers’ Equation 
Another way to treat the nonlinear term in equation (l) is to difference 
the conservation form of Burgers' equation, i.e.. 


1 2 

u. + — (u - u) » vu 

t 2 X XX 


(9) 


The finite difference form of equation (9) is given by equation (2) where now 
' 2 ^'^i+1 ^i-1 with Khosla-Rubin differencing. 

Calculations were made in which a^ was lagged one time step with no 
multiple iteration performed at each time step. The steady-state convergence 
criterion was again e ■ lO"^^. For v = 1/24, converged solutions were 
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obtained for all values of C used (0.9i 10. 100. and lOOO). The solu*> 

max 

tlons are shown in figure For v ■ 1/U6 (R ■ converged solutions 

c imsoc 

were obtained for C • 1000. When v was reduced to 1/96 (R ■ 9*6), 

max c,max 

solutions were obtained for C *0.9* but could not be obtained for 

max 

C > 10. The steady-state solutions for v • 1/U6 and 1/96 are also shown 
max 

on figure 4. No other iteration schemes were tested with the Khosla-Rubin 
method . 

Cases were also run in which a central difference was used for the con- 
vective term. (The procedure of ref. 1 wis used with ^ ^'^i+1 * 1 

lagged one time step and without multiple iteration at each time step.) The 

coefficient matrix is not diagonally dominant when R >2 xinless 

c .meoc 

C <1; however, converged solutions were obtained for v ■ 1/24 (R ■2.4) 
max ~ c ,max 

and for V ■ 1/48 (R =4.8) with values of C up to 1000. Unlike the 

c ,max max 

Khosla-Rubin calculation for v ■ 1/96, the solution converged for = 10, 

although It diverged when C was 100. 

max 

Thus, for Burgers* equation with small values of V, the use of conserva- 
tion form is appropriate. Since the solutions in such cases are steep waves, 
this result is not unexpected. 


CONCLUTING REMARKS 

In previously reported implicit calculations of the nonlinear Burgers' 
equation (in nonconservative form) using central finite differences for 
spatial derivatives, solutions could not be obtained in many (but not all) 
cases when the cell Reynolds nvanber was greater than 2. This failure was 
attributed to the loss of diagonal dominance in the coefficient matrix which 
occurs when Courant numbers greater than one are used for calculations in which 
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the maximum cell Reynolds number exceeds 2. The Khosla-Rubin modificatiou, 
which gives an unconditionally diagonally dcmiinant coefficient matrix, was 
expected to give stable solutions for these cases. Numerical tests, tx>wever, 
have shown that this method is unstable for cases in iidiich the maximum cell 
Reynolds number is greater than 2; i.e., for cases in which **wiggles" are 
inherent in the solution to the difference equation, when the Courant number 
is greater than 1. VHien the Coxirant number was less than 1, stable solutions 
could be obtained for scxne values of cell Reynolds number, but the maximum 
allowable time step was found to be dependent upon th'^ iteration technique 
used. For the linear IVorgers* equation, on the other h^nd, no instability 
was found even when "wigglec” were present. These r«sult'> suggest that the 
failures may be the result of a nonlinear instability. 

In the previous maaerical tests with central differencing, n < cases were 
found in which converged solutions could not be obtained when the maximum 
cell Reynolds number was 2 or less. However, when the maximum cell Reynolds 
number was greater than 2, solutions could not be obtained in scane cases 
even when the time step was small enough to give a diagonally dominant 
coefficient matrix. Here again the failure may have been the result of a 
nonlinear instability. 

This particular instability can be overcome for many cases by using the 
conservative form of Burgers' equation for both the central difference and 
ICiiosla-Bubin metb<^ds. The present study did not, however, attempt to 
determine whether upper limits on the maximum cell Reynolds number and Courant 
number exist for the conservative form. 


11 


REFERENCES 


1. Hirsh t Richard S.» and Rudy, David H. : The Role of Diagonal Dominance and 

Cell Reynolds Number in Implicit Difference Methods for Fluid Mechanics 
Problems, Journal of Ctaanutational Physics , v. l6, no. 3» November 19T*<, 
pp. 30 U 3 IO. 

2. Khosla, P. K.; and Rubin, S. 0.: A Dia^snally D«sinant Second*^rder 

Accurate Implicit Scheme. Computer euid Fluids , v. 2, I 97 I*, 

pp. 207 - 209 . 

3 . Keller, Herbert B.; and Cebeci, Tuncer^ Accurate Numerical Methods for 

Boundary-Layer Flows. II 5 Two-Dimensional Turbulent Flowt , AIAA 
Journal , v. 10, no. 9» September 1972, pp. 1193-1199. 


.4 



TABLE I. SUMMARY OF COMPUTED RESULTS 
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